﻿using System;
using System.Collections.Generic;
using Allegro.Integrator;

namespace Allegro.Integrator.RungeKutta
{
    /// <summary>
    /// This is a stage 9, order (5,6) Runge-Kutta pair published by Jim Verner.
    /// Documentation of this (5,6) Runge-Kutta pair:
    /// http://www.math.sfu.ca/~jverner/FSAL56MostEfficientOptOnWebSite
    /// </summary>
    public class RKTableau_56_Verner01 : ButcherTableau
    {
        public RKTableau_56_Verner01() : base(9)
        {
            // Initialize the vectors
            c[0] = 0;
            c[1] = 3.0 / 50.0;
            c[2] = 1439.0 / 15000.0;
            c[3] = 1439.0 / 10000.0;
            c[4] = 4973.0 / 10000.0;
            c[5] = 389.0 / 400.0;
            c[6] = 1999.0 / 2000.0;
            c[7] = 1.0;
            c[8] = 1.0;

            a[1, 0] = 3.0 / 50.0;

            a[2, 0] = 519479.0 / 27000000.0;
            a[2, 1] = 2070721.0 / 27000000.0;

            a[3, 0] = 1439.0 / 40000.0;
            a[3, 1] = 0.0;
            a[3, 2] = 4317.0 / 40000.0;

            a[4, 0] = 109225017611.0 / 82828840000.0;
            a[4, 1] = 0.0;
            a[4, 2] = -417627820623.0 / 82828840000.0;
            a[4, 3] = 43699198143.0 / 10353605000.0;

            a[5, 0] = -8036815292643907349452552172369.0 / 191934985946683241245914401600.0;
            a[5, 1] = 0.0;
            a[5, 2] = 246134619571490020064824665.0 / 1543816496655405117602368.0;
            a[5, 3] = -13880495956885686234074067279.0 / 113663489566254201783474344.0;
            a[5, 4] = 755005057777788994734129.0 / 136485922925633667082436.0;

            a[6, 0] = -1663299841566102097180506666498880934230261.0 / 30558424506156170307020957791311384232000.0;
            a[6, 1] = 0.0;
            a[6, 2] = 130838124195285491799043628811093033.0 / 631862949514135618861563657970240.0;
            a[6, 3] = -3287100453856023634160618787153901962873.0 / 20724314915376755629135711026851409200.0;
            a[6, 4] = 2771826790140332140865242520369241.0 / 396438716042723436917079980147600.0;
            a[6, 5] = -1799166916139193.0 / 96743806114007800.0;

            a[7, 0] = -832144750039369683895428386437986853923637763.0 / 15222974550069600748763651844667619945204887.0;
            a[7, 1] = 0.0;
            a[7, 2] = 818622075710363565982285196611368750.0 / 3936576237903728151856072395343129.0;
            a[7, 3] = -9818985165491658464841194581385463434793741875.0 / 61642597962658994069869370923196463581866011.0;
            a[7, 4] = 31796692141848558720425711042548134769375.0 / 4530254033500045975557858016006308628092.0;
            a[7, 5] = -14064542118843830075.0 / 766928748264306853644.0;
            a[7, 6] = -1424670304836288125.0 / 2782839104764768088217.0;

            a[8, 0] = 382735282417.0 / 11129397249634.0;
            a[8, 1] = 0.0;
            a[8, 2] = 0.0;
            a[8, 3] = 5535620703125000.0 / 21434089949505429.0;
            a[8, 4] = 13867056347656250.0 / 32943296570459319.0;
            a[8, 5] = 626271188750.0 / 142160006043.0;
            a[8, 6] = -51160788125000.0 / 289890548217.0;
            a[8, 7] = 163193540017.0 / 946795234.0;

            // High order weights
            b[0] = 382735282417.0 / 11129397249634.0;
            b[1] = 0.0;
            b[2] = 0.0;
            b[3] = 5535620703125000.0 / 21434089949505429.0;
            b[4] = 13867056347656250.0 / 32943296570459319.0;
            b[5] = 626271188750.0 / 142160006043.0;
            b[6] = -51160788125000.0 / 289890548217.0;
            b[7] = 163193540017.0 / 946795234.0;
            b[8] = 0.0;

            // Low order weights
            bh[0] = 124310637869885675646798613.0 / 2890072468789466426596827670.0;
            bh[1] = 0.0;
            bh[2] = 0.0;
            bh[3] = 265863151737164990361330921875.0 / 1113197271463372303940319369579.0;
            bh[4] = 3075493557174030806536302953125.0 / 6843749922042323876546949699876.0;
            bh[5] = 67798000008733879813263055.0 / 29532792147666737550036372.0;
            bh[6] = -1099436585155390846238326375.0 / 15055706496446408859196167.0;
            bh[7] = 26171252653086373181571802.0 / 368794478890732346033505.0;
            bh[8] = 1.0 / 30.0;

            //bh[0] = 273361583.0 / 5567482366.0;
            //bh[1] = 0.0;
            //bh[2] = 0.0;
            //bh[3] = 1964687500000.0 / 8727630165387.0;
            //bh[4] = 596054687500.0 / 1269637976277.0;
            //bh[5] = 12740367500.0 / 15795556227.0;
            //bh[6] = 0.0;
            //bh[7] = -4462730789736252634813752317.0 / 7350663039626676022821734166.0;
            //bh[8] = 441454562788983500.0 / 7763730504400359099.0;
        }
    }
}
